###############################################################################
## AUTHOR: ALAN YAN
## DATE LAST UPDATED: 09/27/2020
## PURPOSE: HET EFFECTS VIA INTERACTION (KEY COMPONENT = EDUCATION CONTINUOUS)
###############################################################################
rm(list = ls())

#### NOTES ####

#(1) create a tidier function for random forest output

#### PACKAGES ####
library(pacman)
p_load(tidyverse,
       estimatr,
       broom,
       haven,
       hrbrthemes,
       cjoint,
       formula.tools, 
       DeclareDesign,
       emmeans,
       lmtest,
       cregg)

#### LOAD DATA ####
df <- read_rds("01-yougov-conjoint/01-data/cleaned/yougov-conjoint-cleaned.rds")

#### ANALYZE DATA (NO CONTROLS) ####
#Do people prefer firms with these features?
work <- lm(work_binary ~ 
                                    r_education_continuous*as.factor(cje_corporate_gov) +
                                    r_education_continuous*as.factor(cje_corporate_resp) +
                                    r_education_continuous*as.factor(cje_firm_size) +
                                    r_education_continuous*as.factor(cje_political_donations) + 
                                    r_education_continuous*as.factor(cje_gender_ownership) + 
                                    r_education_continuous*as.factor(cje_healthcare) + 
                                    r_education_continuous*as.factor(cje_hours) + 
                                    r_education_continuous*as.factor(cje_training) + 
                                    r_education_continuous*as.factor(cje_location) + 
                                    r_education_continuous*as.factor(cje_sick_leave) + 
                                    r_education_continuous*as.factor(cje_parental_leave) + 
                                    r_education_continuous*as.factor(cje_race_ownership) + 
                                    r_education_continuous*as.factor(cje_retirement) + 
                                    r_education_continuous*as.factor(cje_income) +
                                    r_education_continuous*as.factor(cje_type_of_work) +
                                    r_education_continuous*as.factor(cje_union) +
                                    r_education_continuous*as.factor(cje_work_from_home) +
                                    r_education_continuous*as.factor(cje_team_work) +
                                    r_education_continuous*as.factor(cje_work_culture) +
                                    r_education_continuous,
                                  data = df
) 

#Do people believe some firms are better at dealing with complaints?
complaints <- lm(complaints_binary ~ 
                                          r_education_continuous*as.factor(cje_corporate_gov) +
                                          r_education_continuous*as.factor(cje_corporate_resp) +
                                          r_education_continuous*as.factor(cje_firm_size) +
                                          r_education_continuous*as.factor(cje_political_donations) + 
                                          r_education_continuous*as.factor(cje_gender_ownership) + 
                                          r_education_continuous*as.factor(cje_healthcare) + 
                                          r_education_continuous*as.factor(cje_hours) + 
                                          r_education_continuous*as.factor(cje_training) + 
                                          r_education_continuous*as.factor(cje_location) + 
                                          r_education_continuous*as.factor(cje_sick_leave) + 
                                          r_education_continuous*as.factor(cje_parental_leave) + 
                                          r_education_continuous*as.factor(cje_race_ownership) + 
                                          r_education_continuous*as.factor(cje_retirement) + 
                                          r_education_continuous*as.factor(cje_income) +
                                          r_education_continuous*as.factor(cje_type_of_work) +
                                          r_education_continuous*as.factor(cje_union) +
                                          r_education_continuous*as.factor(cje_work_from_home) +
                                          r_education_continuous*as.factor(cje_team_work) +
                                          r_education_continuous*as.factor(cje_work_culture) +
                                          r_education_continuous,
                                        data = df
) 

#Do people believe they'll have more responsibilties at some firms?
responsibility <- lm(responsibility_binary ~ 
                                              r_education_continuous*as.factor(cje_corporate_gov) +
                                              r_education_continuous*as.factor(cje_corporate_resp) +
                                              r_education_continuous*as.factor(cje_firm_size) +
                                              r_education_continuous*as.factor(cje_political_donations) + 
                                              r_education_continuous*as.factor(cje_gender_ownership) + 
                                              r_education_continuous*as.factor(cje_healthcare) + 
                                              r_education_continuous*as.factor(cje_hours) + 
                                              r_education_continuous*as.factor(cje_training) + 
                                              r_education_continuous*as.factor(cje_location) + 
                                              r_education_continuous*as.factor(cje_sick_leave) + 
                                              r_education_continuous*as.factor(cje_parental_leave) + 
                                              r_education_continuous*as.factor(cje_race_ownership) + 
                                              r_education_continuous*as.factor(cje_retirement) + 
                                              r_education_continuous*as.factor(cje_income) +
                                              r_education_continuous*as.factor(cje_type_of_work) +
                                              r_education_continuous*as.factor(cje_union) +
                                              r_education_continuous*as.factor(cje_work_from_home) +
                                              r_education_continuous*as.factor(cje_team_work) +
                                              r_education_continuous*as.factor(cje_work_culture) +
                                              r_education_continuous,
                                            data = df
)

#Do people believe they'll have more power at some firms?
power <- lm(power_binary ~ 
              r_education_continuous*as.factor(cje_corporate_gov) +
                                     r_education_continuous*as.factor(cje_corporate_resp) +
                                     r_education_continuous*as.factor(cje_firm_size) +
                                     r_education_continuous*as.factor(cje_political_donations) + 
                                     r_education_continuous*as.factor(cje_gender_ownership) + 
                                     r_education_continuous*as.factor(cje_healthcare) + 
                                     r_education_continuous*as.factor(cje_hours) + 
                                     r_education_continuous*as.factor(cje_training) + 
                                     r_education_continuous*as.factor(cje_location) + 
                                     r_education_continuous*as.factor(cje_sick_leave) + 
                                     r_education_continuous*as.factor(cje_parental_leave) + 
                                     r_education_continuous*as.factor(cje_race_ownership) + 
                                     r_education_continuous*as.factor(cje_retirement) + 
                                     r_education_continuous*as.factor(cje_income) +
                                     r_education_continuous*as.factor(cje_type_of_work) +
                                     r_education_continuous*as.factor(cje_union) +
                                     r_education_continuous*as.factor(cje_work_from_home) +
                                     r_education_continuous*as.factor(cje_team_work) +
                                     r_education_continuous*as.factor(cje_work_culture) +
                                     r_education_continuous,
                                   data = df
)

#### tables S64 to S87 ####
dep_var_labels <- c("Prefer to work",
                    "More power",
                    "More responsibilities",
                    "Better handle complaints")
covariate_labels <- names(work$coefficients)[-1] %>%
                        stringr::str_replace("r_education_continuous", "Education") %>%
                        stringr::str_replace(":", " X ") %>%
                        stringr::str_remove_all("as.factor\\(cje_corporate_gov\\)|as.factor\\(cje_corporate_resp\\)|as.factor\\(cje_firm_size\\)|as.factor\\(cje_political_donations\\)|as.factor\\(cje_gender_ownership\\)|as.factor\\(cje_healthcare\\)|as.factor\\(cje_hours\\)|as.factor\\(cje_training\\)|as.factor\\(cje_location\\)|as.factor\\(cje_sick_leave\\)|as.factor\\(cje_parental_leave\\)|as.factor\\(cje_race_ownership\\)|as.factor\\(cje_retirement\\)|as.factor\\(cje_income\\)|as.factor\\(cje_type_of_work\\)|as.factor\\(cje_union\\)|as.factor\\(cje_work_from_home\\)|as.factor\\(cje_team_work\\)|as.factor\\(cje_work_culture\\)")

stargazer(
  work,
  power,
  dep.var.labels = dep_var_labels[1:2],
  covariate.labels = covariate_labels,
  se = list(
    starprep(work, stat = "std.error", clusters = df$rid, se_type = "stata")[[1]],
    starprep(power, stat = "std.error", clusters = df$rid, se_type = "stata")[[1]]
  ),
  p = list(
    starprep(work, stat = "p.value", clusters = df$rid, se_type = "stata")[[1]],
    starprep(power, stat = "p.value", clusters = df$rid, se_type = "stata")[[1]]
  ),
  title = "Replicating Figure A7 work and power results",
  label = "tab:education-interaction-main-effects-work-power"
)

stargazer(
  responsibility,
  complaints,
  dep.var.labels = dep_var_labels[3:4],
  covariate.labels = covariate_labels,
  se = list(
    starprep(responsibility, stat = "std.error", clusters = df$rid, se_type = "stata")[[1]],
    starprep(complaints, stat = "std.error", clusters = df$rid, se_type = "stata")[[1]]
  ),
  p = list(
    starprep(responsibility, stat = "p.value", clusters = df$rid, se_type = "stata")[[1]],
    starprep(complaints, stat = "p.value", clusters = df$rid, se_type = "stata")[[1]]
  ),
  title = "Replicating Figure A7 responsibility and complaint results",
  label = "tab:education-interaction-main-effects-responsibility-complaints"
)
